• Wstęp
  • Wczytanie bibliotek
  • Wczytanie danych.
  • Czyszczenie zbioru.
  • Podsumowanie zbioru.
  • Wpływ wartości wybranych atrybutów na zgon pacjenta.
  • Pacjenci opuszczający szpital w danym dniu.
  • Przygotowanie klasyfikatora przewidującego czy dany pacjent przeżyje.
  • Uczenie klasyfikatora.

Wstęp

Projekt zakładał analizę danych pacjentów chorych na COVID-19. Przetworzenie danych ze zbioru umożliwiło sporządzenie raportu zawierającego wykresy takie jak rozkład wyzdrowień i śmierci. Podjęto także próbę stworzenia klasyfikatora, który będzie potrafił wskazać czy dany pacjent przeżyje.

Wczytanie bibliotek

library("readxl")
library("dplyr")
library(skimr)
library(plotly)
library(ggplot2)
library(DT)
library(hablar)
library(reshape2)
library(tidyr)
library(caret)
set.seed(23)

Wczytanie danych.

data = data.frame(read_excel("wuhan_blood_sample_data_Jan_Feb_2020.xlsx"))

Czyszczenie zbioru.

Data próbki została przycięte do formatu “yyyy-MM-dd”.Następnie dane zostały zgrupowane pod dacie przyjęcia oraz dacie opuszenia szpitala, z powodu wyzdrowienia lub zgonu. Co umożliwiło przypisanie brakujących PATIENT_ID.

df <- mutate(data, RE_DATE = as.Date(substr(RE_DATE,1,10), "%Y-%m-%d"))
df <- df %>% group_by(Admission.time, Discharge.time) %>% mutate(PATIENT_ID = max(PATIENT_ID, na.rm = TRUE)) %>% ungroup()
df <- df %>% group_by(PATIENT_ID, RE_DATE) %>% 
  mutate(Admission.time = max(Admission.time, na.rm = TRUE)) %>%
  mutate(Discharge.time = max(Discharge.time, na.rm = TRUE))%>% ungroup()

Wszystkie próbki z danego dnia dla danego pacjenta zostały zgrupowane.

df <- df %>% group_by(PATIENT_ID, RE_DATE, Admission.time, Discharge.time, age, outcome, gender) %>% summarise_all(sum_) %>% ungroup()

Dane liczbowe odpowiadające płci zostały zamienione na dane tekstowe.

df <- df %>% mutate(gender = factor(gender, labels = c('male', 'female')))
df <- df %>% mutate(outcome = factor(outcome, labels = c('survived', 'died')))

Kolumny i wiersze zawierające niekompletne dane zostały usunięte.

df$X2019.nCoV.nucleic.acid.detection = NULL
df <- subset(df, !is.na(RE_DATE))

Tabela przedstawia zbiór po przetworzeniu.

datatable(df,options = list(scrollX = TRUE), filter = 'top')

Podsumowanie zbioru.

my_skim <- skim_with(character = sfl(), append = FALSE, base = sfl(n_missing))
datatable(my_skim(df), options = list(scrollX = TRUE))

Wpływ wartości wybranych atrybutów na zgon pacjenta.

Wykresy przedstawiają wpływ wartości dehydrogenazy mleczanowej (LDH), białka C-reaktywnego (hs-CRP), oraz liczby limfocytów.

last_data <- df
last_data <- last_data %>% group_by(PATIENT_ID) %>%
  summarise(outcome = last(outcome), LDH = last(Lactate.dehydrogenase), Limfocyty = last(Lactate.dehydrogenase), 'hs-CRP' = last(High.sensitivity.C.reactive.protein)) %>% ungroup()

last_data[ ,c('PATIENT_ID')] <- list(NULL)
datatable(last_data)
last_plot_data = melt(last_data, id=c("outcome"))
graph <- ggplot(last_plot_data, aes(x=value, y=outcome)) + geom_point() + facet_wrap(~variable)
ggplotly(graph)
0100020003000surviveddied01000200030000100020003000
valueoutcomeLDHLimfocytyhs-CRP

Pacjenci opuszczający szpital w danym dniu.

Wykres przedstawia dane dla każdego dnia. Zsumowana została liczba pacjentów opuszczających szpital z podziałem na zgony oraz wyzdrowienia.

by_date <- mutate(df, Discharge.time = substr(Discharge.time,1,10))
by_date <- by_date %>% group_by(Discharge.time, PATIENT_ID) %>% mutate(outcome = first(outcome == "died")) %>%
  summarise(outcome = first(outcome)) %>% group_by(Discharge.time) %>%
  summarise(Pacjenci = n_distinct(PATIENT_ID), Zgony = sum(outcome == TRUE), Wyzdrowienia = sum(outcome == FALSE))

plot_data = melt(by_date, id=c("Discharge.time"))

p <- ggplot(plot_data) + 
  geom_point(aes(x = Discharge.time, y = value, color = variable,
                 text = paste(' Dzień: ', Discharge.time,'<br>','Liczba: ',value,'<br>','Typ: ',variable))) + 
  scale_colour_manual(values=c("black","red","green")) +
  labs(color = "Legenda", x = "Dzień", y = "Liczba") +
  theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1))

g <- ggplotly(p, tooltip = "text") %>% layout(legend = list(x = 0.05, y = 0.95))
g
2020-01-232020-01-292020-01-302020-01-312020-02-012020-02-022020-02-032020-02-042020-02-052020-02-062020-02-072020-02-082020-02-092020-02-102020-02-112020-02-122020-02-132020-02-142020-02-152020-02-162020-02-172020-02-182020-02-192020-02-202020-02-262020-03-012020-03-022020-03-032020-03-0401020304050
PacjenciZgonyWyzdrowieniaDzieńLiczbaLegenda

Przygotowanie klasyfikatora przewidującego czy dany pacjent przeżyje.

Przed przystąpienie do uczenia klasyfikatora dane zostały odpowiednio przygotowane. Ze zbioru zostali usunięci pacjenci ze zbyt mała liczbą przeprowadzonych badan oraz badańia wykonane na zbyt małej liczbie pacjentów.

raw_train_data <- df
raw_train_data[ ,c('Admission.time', 'Discharge.time')] <- list(NULL)
raw_train_data <- raw_train_data[which(rowMeans(!is.na(raw_train_data)) > 0.7), which(colMeans(!is.na(raw_train_data)) > 0.32)]

Po wstępnym przetworzeniu zbioru, usunięto wszystkich pacjentów, którzy nie mieli kompletu pozostałych badań.

raw_train_data <- raw_train_data %>% drop_na()
patient_data <- raw_train_data %>%  group_by(PATIENT_ID) %>% summarise_at(vars(outcome, gender), last) %>% ungroup()
raw_train_data <- raw_train_data %>%  group_by(PATIENT_ID) %>% summarise_if(is.numeric, max) %>% ungroup()

Tabela przedstawia ostateczny zbiór, który zostanie wykorzystany do uczenia klasyfikatora. Zawiera on 217 rekordóW.

raw_train_data <- merge(x = patient_data, y = raw_train_data, by ="PATIENT_ID", all = TRUE)
datatable(raw_train_data,options = list(scrollX = TRUE))

Dane treningowe stanowią 75% zbioru. Pozostałe 25% zostało podzielone pomiędzy dane testowe (40%) i walidujące (60%).

inTraining <- createDataPartition(y = raw_train_data$outcome, p = .75, list = FALSE)
training <- raw_train_data[inTraining,]
testing_val  <- raw_train_data[-inTraining,]

testing_val_index <- createDataPartition(y = testing_val$outcome, p = .6, list = FALSE)
testing <- testing_val[-testing_val_index,]
val <- testing_val[testing_val_index,]

datatable(training, options = list(scrollX = TRUE))
datatable(testing, options = list(scrollX = TRUE))
datatable(val, options = list(scrollX = TRUE))

Uczenie klasyfikatora.

Dane zostały wykorzystane do uczenia klasyfikator random forest.

ctrl <- trainControl(method = "repeatedcv", number = 2,repeats = 5)
fit <- train(outcome ~ .,
             data = training,
             method = "rf",
             trControl = ctrl,
             ntree = 10)

fit
## Random Forest 
## 
## 163 samples
##  58 predictor
##   2 classes: 'survived', 'died' 
## 
## No pre-processing
## Resampling: Cross-Validated (2 fold, repeated 5 times) 
## Summary of sample sizes: 81, 82, 81, 82, 82, 81, ... 
## Resampling results across tuning parameters:
## 
##   mtry  Accuracy   Kappa    
##    2    0.9164408  0.8305415
##   30    0.9839958  0.9676142
##   58    0.9950918  0.9901007
## 
## Accuracy was used to select the optimal model using the largest value.
## The final value used for the model was mtry = 58.
rfGrid <- expand.grid(mtry = 10:30)
gridCtrl <- trainControl(method = "repeatedcv", summaryFunction = twoClassSummary, classProbs = TRUE, number = 2, repeats = 5)
fitTune <- train(outcome ~ .,
             data = training,
             method = "rf",
             metric = "ROC",
             preProc = c("center", "scale"),
             trControl = gridCtrl,
             tuneGrid = rfGrid,
             ntree = 30)
fitTune
## Random Forest 
## 
## 163 samples
##  58 predictor
##   2 classes: 'survived', 'died' 
## 
## Pre-processing: centered (58), scaled (58) 
## Resampling: Cross-Validated (2 fold, repeated 5 times) 
## Summary of sample sizes: 81, 82, 82, 81, 81, 82, ... 
## Resampling results across tuning parameters:
## 
##   mtry  ROC        Sens       Spec     
##   10    0.9933192  0.9666667  0.9669670
##   11    0.9986728  0.9777778  0.9835586
##   12    0.9976852  0.9688889  0.9807057
##   13    0.9978095  0.9755556  0.9863363
##   14    0.9971330  0.9688889  0.9725976
##   15    0.9971296  0.9755556  0.9694444
##   16    0.9979630  0.9844444  0.9780030
##   17    0.9983025  0.9755556  0.9888889
##   18    0.9981790  0.9666667  1.0000000
##   19    0.9982716  0.9777778  0.9917417
##   20    0.9991358  0.9822222  0.9945195
##   21    0.9984259  0.9733333  0.9833333
##   22    0.9989506  0.9888889  0.9889640
##   23    0.9983959  0.9844444  0.9972973
##   24    0.9988889  0.9822222  0.9889640
##   25    0.9992009  0.9866667  1.0000000
##   26    0.9990741  0.9844444  0.9972222
##   27    0.9999074  0.9866667  1.0000000
##   28    0.9991975  0.9844444  1.0000000
##   29    0.9996914  0.9888889  0.9972973
##   30    0.9992593  0.9888889  1.0000000
## 
## ROC was used to select the optimal model using the largest value.
## The final value used for the model was mtry = 27.

Uzyskano bardzo wysoki wskaźnik “accuracy”. Nauczony klasyfikator poprawnie sklasyfikował wszystkie próbki. Niestety liczba danych wykorzystana w procesie jest za mała. Nie można mieć pewności, że klasyfikator będzie działać z tak dobrą skutecznością na znacznie większym zbiorze.

predict_result <- predict(fitTune, newdata = val)
confusionMatrix(data = predict_result, val$outcome)
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction survived died
##   survived       18    0
##   died            0   15
##                                      
##                Accuracy : 1          
##                  95% CI : (0.8942, 1)
##     No Information Rate : 0.5455     
##     P-Value [Acc > NIR] : 2.056e-09  
##                                      
##                   Kappa : 1          
##                                      
##  Mcnemar's Test P-Value : NA         
##                                      
##             Sensitivity : 1.0000     
##             Specificity : 1.0000     
##          Pos Pred Value : 1.0000     
##          Neg Pred Value : 1.0000     
##              Prevalence : 0.5455     
##          Detection Rate : 0.5455     
##    Detection Prevalence : 0.5455     
##       Balanced Accuracy : 1.0000     
##                                      
##        'Positive' Class : survived   
##